Análisis de conjuntos frecuentes

Una de las tareas más antiguas de la minería de datos es la búsqueda de conjuntos frecuentes en canastas, o un análisis derivado que se llama análisis de reglas de asociación.

Originalmente, pensamos que tenemos una colección grande de tickets de un supermercado. Nos interesa encontrar subconjuntos de artículos (por ejemplo, pan y leche) que ocurren frecuentemente en esos tickets. La idea es que si tenemos estos subconjuntos frecuentes, entonces podemos entender mejor el tipo de compras que hacen los clientes, diseñar mejor promociones y entender potenciales efectos cruzados, reordenar los estantes del supermercado, etc.

En general, los conjuntos frecuentes indican asociaciones (y cuantificaciones de la asociación) entre artículos que hay que tomar en cuenta al momento de tomar decisiones. Esto normalmente se llama análisis de market basket.

El análisis de subconjuntos frecuentes puede ser utilizado para otros propósitos, como veremos más adelante.

Datos de canastas

Consideremos el siguiente ejemplo chico del paquete arules. Contiene unas 10 mil canastas observadas en de un supermercado durante un mes, agregadas a 169 categorías.

data(Groceries) # del paquete arules
class(Groceries)

## [1] "transactions"
## attr(,"package")
## [1] "arules"

head(Groceries)

## transactions in sparse format with
##  6 transactions (rows) and
##  169 items (columns)

lista_mb <- as(Groceries, 'list')

Estos tres canastas (tickets) de ejemplo:

lista_mb[[2]]
## [1] "tropical fruit" "yogurt"         "coffee"
lista_mb[[52]]
## [1] "canned beer"
lista_mb[[3943]]
## [1] "sausage"        "UHT-milk"       "flour"          "flower (seeds)"

Podemos calcular la distribución del número de artículos por canasta:

sprintf("Número de canastas: %s", length(lista_mb))

## [1] "Número de canastas: 9835"

num_items <- sapply(lista_mb, length)
sprintf("Promedio de artículos por canasta: %.3f", mean(num_items))

## [1] "Promedio de artículos por canasta: 4.409"

qplot(num_items, binwidth=1)

Podemos hacer una tabla con las canastas y examinar los artículos más frecuentes:

canastas_nest <- data_frame(canasta_id = 1:length(lista_mb),
                      articulos = lista_mb)
head(canastas_nest)

## # A tibble: 6 x 2
##   canasta_id articulos
##        <int> <list>   
## 1          1 <chr [4]>
## 2          2 <chr [3]>
## 3          3 <chr [1]>
## 4          4 <chr [4]>
## 5          5 <chr [4]>
## 6          6 <chr [5]>

canastas <- canastas_nest %>% unnest
canastas

## # A tibble: 43,367 x 2
##    canasta_id articulos          
##         <int> <chr>              
##  1          1 citrus fruit       
##  2          1 semi-finished bread
##  3          1 margarine          
##  4          1 ready soups        
##  5          2 tropical fruit     
##  6          2 yogurt             
##  7          2 coffee             
##  8          3 whole milk         
##  9          4 pip fruit          
## 10          4 yogurt             
## # ... with 43,357 more rows

num_canastas <- nrow(canastas_nest)
num_canastas

## [1] 9835

articulos_frec <- canastas %>% group_by(articulos) %>%
                  summarise(n  = n()) %>%
                  mutate(prop = n / num_canastas) %>%
                  arrange(desc(n))
DT::datatable(articulos_frec %>%
  mutate_if(is.numeric, funs(round(., 3))))

Un primer análisis que podríamos considerar es el de canastas completas que ocurren frecuentemente. ¿Qué tan útil crees que puede ser este análisis?

colapsar_canasta <- function(x, sep = '-'){
  # convierte cada canasta a una cadena
  x %>% as.character %>% sort %>% paste(collapse = '-')
}
canastas_conteo <- canastas_nest %>%
                rowwise() %>%
                mutate(canasta_str = colapsar_canasta(articulos)) %>%
                group_by(canasta_str) %>%
                summarise(n = n()) %>%
                mutate(prop = round(n /num_canastas, 5)) %>%
                arrange(desc(n))
nrow(canastas_conteo)

## [1] 7011

Y aquí vemos las canastas más frecuentes:

DT::datatable(canastas_conteo %>% head(n = 100) %>%
    mutate_if(is.numeric, funs(round(., 4))))

Hay algunas canastas (principalmente canastas que contienen un solo artículo) que aparecen con frecuencia considerable (alrededor de 1% o 2%), pero las canastas están bastante dispersas en el espacio de posibles canastas (que es gigantesco: ¿puedes calcularlo?). Debido a esta dispersión (las canastas se distribuyen de manera rala en el espacio de posibles canastas), este análisis es de utilidad limitada.

Datos de canastas

El tamaño de las canastas normalmente es muy chico (por ejemplo 1 a 30 artículos distintos).
El número de artículos típicamente no es muy grande (de cientos a cientos de miles, por ejemplo).
El número de canastas puede ser mucho mayor (en algunos casos miles de millones) y quizá no pueden leerse completas en memoria.
El número de canastas distintas es alto, y hay pocas canastas frecuentes.

El último inciso señala que encontrar canastas frecuentes no será muy informativo. En lugar de eso buscamos conjuntos de artículos (que podríamos llamar subcanastas) que forman parte de muchas canastas.

Conjuntos frecuentes

Un enfoque simple y escalable para analizar estas canastas es el de los conjuntos frecuentes (frequent itemsets).

Conjuntos frecuentes

Consideramos un conjunto de artículos \(I = \{s_1,s_2,\ldots, s_k\}\). El soporte de \(I\) lo definimos como la proporción de canastas que contienen (al menos) estos artículos: \[P(I) = \frac{n(I)}{n},\] donde \(n(I)\) es el número de canastas que contienen todos los artículos de \(I\), y \(n\) es el número total de canastas.

Sea \(s\in (0,1)\). Para este valor fijo \(s\), decimos que un conjunto de artículos \(I\) es un conjunto frecuente cuando \(P(I)\geq s\).

Observación: Si hay \(m\) artículos, entonces el número de posibles itemsets es de \(2^m\). El número total de itemsets de tamaño \(k\) es \(\binom{m}{k}\). Sin embargo, esperamos que los itemsets frecuentes sean raros, de forma que es factible calcularlos (siempre y cuando el soporte mínimo sea suficientemente grande).

Ejercicio

Considera las canastas {1,2,3}, {1,2}, {2,4}, {2,3}. ¿Cuáles son los itemsets frecuentes de soporte > 0.4?

lista <- list()
lista[[1]] <- c(1,2,3)
lista[[2]] <- c(1,2)
lista[[3]] <- c(2,4)
lista[[4]] <- c(2,3)

lista

## [[1]]
## [1] 1 2 3
## 
## [[2]]
## [1] 1 2
## 
## [[3]]
## [1] 2 4
## 
## [[4]]
## [1] 2 3

pars <- list(supp = 0.4, target='frequent itemsets')
ap <- apriori(lista, parameter = pars)

## Apriori
## 
## Parameter specification:
##  confidence minval smax arem  aval originalSupport maxtime support minlen
##          NA    0.1    1 none FALSE            TRUE       5     0.4      1
##  maxlen            target   ext
##      10 frequent itemsets FALSE
## 
## Algorithmic control:
##  filter tree heap memopt load sort verbose
##     0.1 TRUE TRUE  FALSE TRUE    2    TRUE
## 
## Absolute minimum support count: 1 
## 
## set item appearances ...[0 item(s)] done [0.00s].
## set transactions ...[4 item(s), 4 transaction(s)] done [0.00s].
## sorting and recoding items ... [3 item(s)] done [0.00s].
## creating transaction tree ... done [0.00s].
## checking subsets of size 1 2 done [0.00s].
## writing ... [5 set(s)] done [0.00s].
## creating S4 object  ... done [0.00s].

inspect(ap)

##     items support count
## [1] {1}   0.5     2    
## [2] {3}   0.5     2    
## [3] {2}   1.0     4    
## [4] {1,2} 0.5     2    
## [5] {2,3} 0.5     2

Ejemplo

Explicamos más adelante la función apriori de arules, pero por lo pronto podemos examinar algunos conjuntos frecuentes de soporte mínimo 0.005 (como hay alrededor de 10000 canastas, esto significa que canastas que aparecieron al menos 50 veces durante el mes):

pars <- list(supp = 0.005, target='frequent itemsets')
ap <- apriori(lista_mb, parameter = pars)

class(ap)

## [1] "itemsets"
## attr(,"package")
## [1] "arules"

head(ap)

## set of 6 itemsets

size(ap)

##    [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##   [35] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##   [69] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
##  [103] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [137] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [171] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [205] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [239] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [273] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [307] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [341] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [375] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [409] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [443] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [477] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [511] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [545] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [579] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [613] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [647] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [681] 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
##  [715] 2 2 2 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [749] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [783] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [817] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [851] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [885] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [919] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [953] 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3 3
##  [987] 3 3 3 4 4 4 4 4 4 4 4 4 4 4 4

Veamos algunos conjuntos frecuentes de tamaño 1:

ap_1 <- subset(ap, size(ap) == 1) 
length(ap_1)
## [1] 120

inspect(sort(ap_1, by='support')[1:10])

##      items              support    count
## [1]  {whole milk}       0.25551601 2513 
## [2]  {other vegetables} 0.19349263 1903 
## [3]  {rolls/buns}       0.18393493 1809 
## [4]  {soda}             0.17437722 1715 
## [5]  {yogurt}           0.13950178 1372 
## [6]  {bottled water}    0.11052364 1087 
## [7]  {root vegetables}  0.10899847 1072 
## [8]  {tropical fruit}   0.10493137 1032 
## [9]  {shopping bags}    0.09852567  969 
## [10] {sausage}          0.09395018  924

Algunas de tamaño 2 y 3:

ap_2 <- subset(ap, size(ap) == 2)
length(ap_2)
## [1] 605

inspect(sort(ap_2, by='support')[1:10])

##      items                              support    count
## [1]  {other vegetables,whole milk}      0.07483477 736  
## [2]  {rolls/buns,whole milk}            0.05663447 557  
## [3]  {whole milk,yogurt}                0.05602440 551  
## [4]  {root vegetables,whole milk}       0.04890696 481  
## [5]  {other vegetables,root vegetables} 0.04738180 466  
## [6]  {other vegetables,yogurt}          0.04341637 427  
## [7]  {other vegetables,rolls/buns}      0.04260295 419  
## [8]  {tropical fruit,whole milk}        0.04229792 416  
## [9]  {soda,whole milk}                  0.04006101 394  
## [10] {rolls/buns,soda}                  0.03833249 377

ap_3 <- subset(ap, size(ap) == 3)
length(ap_3)
## [1] 264

inspect(sort(ap_3, by='support')[1:5])

##     items                                         support    count
## [1] {other vegetables,root vegetables,whole milk} 0.02318251 228  
## [2] {other vegetables,whole milk,yogurt}          0.02226741 219  
## [3] {other vegetables,rolls/buns,whole milk}      0.01789527 176  
## [4] {other vegetables,tropical fruit,whole milk}  0.01708185 168  
## [5] {rolls/buns,whole milk,yogurt}                0.01555669 153

También podemos ver qué itemsets incluyen algún producto particular, por ejemplo

ap_berries <- subset(ap, items %pin% 'berries')
length(ap_berries)
## [1] 10

inspect(sort(ap_berries, by ='support')[1:5])

##     items                        support     count
## [1] {berries}                    0.033248602 327  
## [2] {berries,whole milk}         0.011794611 116  
## [3] {berries,yogurt}             0.010574479 104  
## [4] {berries,other vegetables}   0.010269446 101  
## [5] {berries,whipped/sour cream} 0.009049314  89

ap_soda <- subset(ap, items %pin% 'soda')
length(ap_soda)
## [1] 90

inspect(sort(ap_soda, by ='support')[1:5])

##     items                   support    count
## [1] {soda}                  0.17437722 1715 
## [2] {soda,whole milk}       0.04006101  394 
## [3] {rolls/buns,soda}       0.03833249  377 
## [4] {other vegetables,soda} 0.03274021  322 
## [5] {bottled water,soda}    0.02897814  285

Monotonicidad de conjuntos frecuentes

Ahora consideramos el problema de encontrar los conjuntos frecuentes.

Como discutimos arriba en las características de los datos de canastas, el número de transacciones puede ser muy grande (hasta miles de millones), las canastas típicamente son chicas, y el número de artículos distintos puede ir de las decenas hasta cientos de miles o algunos millones. Los algoritmos que se utilizan están diseñado para tratar con estas características. En particular, suponemos que

La lista de transacciones es muy grande, y no puede leerse completa en memoria,
Sin embargo, para una sola canasta, es posible calcular de manera relativamente rápida todos los subconjuntos de tamaño \(k\) (para \(k=1,2,3,4\), por ejemplo). Por ejemplo, si una canasta tiene 10 artículos, hay \(\binom{10}{3}\) canastas de tamaño 3, \(\binom{10}{4} = 210\). Calcular estos subconjuntos es relativamente rápido comparado con leer de disco una transacción.
Finalmente, suponemos que el número de itemsets frecuentes es relativamente chico (lo cual depende de que escojamos un soporte suficientemente alto).

El principio básico que hace posible hacer los conteos de itemsets frecuentes es el siguiente:

Monotonicidad de itemsets Sea \(s\) un nivel de soporte mínimo fijo.

Si un itemset \(I\) es frecuente, entonces todos sus subconjuntos son itemsets frecuentes.
Equivalentemente, si algún subconjunto de un itemset no es frecuente, entonces el itemset no puede ser frecuente. Así que a priori, no es necesario examinar o contar itemsets que contienen al menos un subconjunto que no sea frecuente.

Este hecho, junto con la selección de un soporte mínimo para los itemsets frecuentes es el que hace que la búsqueda y conteo de itemsets frecuentes sea un problema factible.

La demostración es como sigue: Sea \(n(I)\) el número de canastas que contiene \(I\), y supongamos que \(n(I)>s\) (\(I\) es un conjunto frecuente). Sea ahora \(J\subset I\). Entonces cualquier canasta que contiene los artículos de \(I\) contiene también los artículos de \(J\) (que son menos), de forma que \(n(J)\geq n(I)\). Como \(n(I)>s\), entonces \(J\) es un conjunto frecuente.

Observación: Este hecho a veces es un poco confuso (es monotonicidad decreciente), como demuestra este ejemplo de Kahneman: Linda is 31 years old, single, outspoken, and very bright. She majored in philosophy. As a student, she was deeply concerned with issues of discrimination and social justice, and also participated in anti-nuclear demonstrations. Which is more probable?

Linda is a bank teller.
Linda is a bank teller and is active in the feminist movement.

Algoritmo a-priori

Para entender cómo aplicamos monotonicidad, consideremos cómo calcularíamos los pares frecuentes.

Primero calculamos los artículos frecuentes (itemsets de tamaño 1), que son los artículos que aparecen en al menos una proporción \(s\) de las canastas.

(Contar candidatos) Esto requiere recorrer el archivo de transacciones y contar todos los artículos.
(Podar) Examinamos los conteos y seleccionamos aquellos artículos que son frecuentes.

Por el principio de monotonicidad, ningún par frecuente puede contener un artículo no frecuente. Así que para calcular pares:

(Contar candidatos) Recorremos el archivo de transacciones. Para cada transacción, solo contamos pares candidatos cuyos dos artículos son artículos frecuentes (del paso anterior)
(Podar) Examinamos los conteos y seleccionamos aquellos pares que son frecuentes.

Nótese que este algoritmo requiere dos pasadas sobre el conjunto de transacciones.

Ejercicio

Aplica este algoritmo para las canastas {1,2,3}, {1,8}, {2,4}, {2,3,6,7,8}, {2,3,8}, {1,3,8}, {1,2,3,5}, {2,3}. (soporte > 0.25)

lista <- list()
lista[[1]] <- c(1,2,3)
lista[[2]] <- c(1,8)
lista[[3]] <- c(2,4)
lista[[4]] <- c(2,3,6,7,8)
lista[[5]] <- c(2,3,8)
lista[[6]] <- c(1,3,8)
lista[[7]] <- c(1,2,3,5)
lista[[8]] <- c(2,3)

lista

## [[1]]
## [1] 1 2 3
## 
## [[2]]
## [1] 1 8
## 
## [[3]]
## [1] 2 4
## 
## [[4]]
## [1] 2 3 6 7 8
## 
## [[5]]
## [1] 2 3 8
## 
## [[6]]
## [1] 1 3 8
## 
## [[7]]
## [1] 1 2 3 5
## 
## [[8]]
## [1] 2 3

pars <- list(supp = 0.25, target='frequent itemsets')
ap <- apriori(lista, parameter = pars)

## Apriori
## 
## Parameter specification:
##  confidence minval smax arem  aval originalSupport maxtime support minlen
##          NA    0.1    1 none FALSE            TRUE       5    0.25      1
##  maxlen            target   ext
##      10 frequent itemsets FALSE
## 
## Algorithmic control:
##  filter tree heap memopt load sort verbose
##     0.1 TRUE TRUE  FALSE TRUE    2    TRUE
## 
## Absolute minimum support count: 2 
## 
## set item appearances ...[0 item(s)] done [0.00s].
## set transactions ...[8 item(s), 8 transaction(s)] done [0.00s].
## sorting and recoding items ... [4 item(s)] done [0.00s].
## creating transaction tree ... done [0.00s].
## checking subsets of size 1 2 3 done [0.00s].
## writing ... [12 set(s)] done [0.00s].
## creating S4 object  ... done [0.00s].

inspect(ap)

##      items   support count
## [1]  {1}     0.500   4    
## [2]  {8}     0.500   4    
## [3]  {2}     0.750   6    
## [4]  {3}     0.750   6    
## [5]  {1,8}   0.250   2    
## [6]  {1,2}   0.250   2    
## [7]  {1,3}   0.375   3    
## [8]  {2,8}   0.250   2    
## [9]  {3,8}   0.375   3    
## [10] {2,3}   0.625   5    
## [11] {1,2,3} 0.250   2    
## [12] {2,3,8} 0.250   2

Observaciones

En este algoritmo, no es necesario leer todas las transacciones a la vez, podemos procesarlas por bloques, por ejemplo.
En la primera pasada del algoritmo (artículos frecuentes), no es un problema mantener el conteo de todos los artículos en memoria.
Si en la segunda pasada no usáramos monotonicidad, tendríamos que mantener conteos de todos los posibles pares (que son del orden \(m^2\), donde \(m\) es el número de artículos). Mantener este conteo en memoria podría ser difícil si el número de artículos es grande. Sin embargo, el número de artículos frecuentes generalmente es considerablemente menor.

Para itemsets de tamaño más grande, el algoritmo original a priori [@Agrawal] es:

Algoritmo a-priori

Sea \(L_1\) el número de itemsets frecuentes de tamaño \(1\).

Para obtener \(L_k\), el número de itemsets frecuentes de tamaño \(k\):

Sea \(C_k\) el conjunto de candidatos, construido a partir de \(L_{k-1}\)
Para cada transacción \(t\),

Calculamos \(C_t\), que son los candidatos en \(C_k\) que están en \(t\).
Agregamos 1 a cada conteo de los candidatos en \(C_t\).

Filtramos los elementos de \(C_k\) que tengan conteo mayor que el soporte definido para obtener \(L_k\)
Seguimos hasta que encontramos que algún \(L_k\) es vacío (no hay itemset frecuente)

Observaciones: Este algoritmo se puede implementar de distintas maneras, por ejemplo:

Hay distintas maneras de generar el conjunto \(C_{k}\) de candidatos. El paper original sugiere (suponiendo que los artículos siempre están ordenados en los itemsets) hacer un join de \(L_{k-1}\) consigo misma. Por ejemplo, para generar los tríos \(C_3\) a partir de \(L_2\) hacemos

SELECT a.item1, a.item2, b.item2
FROM L2 AS a, L2 AS b
WHERE a.item1 = a.item1, a.item2 < b.item2

donde es crucial que los itemsets estén ordenados (por índice o lexicográficamente).

Hay también distintas maneras de calcular \(C_t\) para cada transacción. El paper original sugiere una estructura de árbol para encontrar los subconjuntos de \(t\) que están en \(C_k\). Ver también [@borgelt].
Más detalles de la implementación de los algoritmos (incluyendo algunos más modernos como FPGrowth, que está implementado en spark) se puede encontrar en [@mmd] y en [@borgelt].

Modelos simples para análisis de canastas

Podemos entender mejor el comportamiento de este análisis con algunos modelos simples para datos de canastas.

En primer lugar, pensamos que los datos están en forma de codificación dummy (aunque no usemos esta representación para los datos reales, podemos considerarlo teóricamente). Una canasta es entonces un renglón de ceros y unos, dependendiendo qué artículos están o no en la canasta:

\[X= (X_1,X_2,\ldots, X_m)\]

donde \(X_i = 1\) si el artículo \(i\) está en la canasta, y \(X_i=0\) si no está.

Podríamos pensar entonces en construir modelos para la conjunta de las canastas

\[P(X_1=x_1,X_2=x_2,\ldots, X_m=x_m)\]

Ejemplo

Por ejemplo, si los items son 1-camisa, 2-pantalones, 3-chamarra, podríamos tener las dos transacciones

\(X = (1,0,0)\), para alguien que solo compró una camisa
\(X = (1,0,1)\), para alguien que solo compró camisa y chamarra

Y podemos inventar una conjunta para todas las canastas, por ejemplo

##   p c ch
## 1 0 0  0
## 2 1 0  0
## 3 0 1  0
## 4 1 1  0
## 5 0 0  1
## 6 1 0  1
## 7 0 1  1
## 8 1 1  1

set.seed(280)
x <-  c(0, rpois(7, 4)^2)
x

## [1]  0  9 16 16 16 16  9  4

probs$prob <- round(x/sum(x),3)
probs

##   p c ch  prob
## 1 0 0  0 0.000
## 2 1 0  0 0.105
## 3 0 1  0 0.186
## 4 1 1  0 0.186
## 5 0 0  1 0.186
## 6 1 0  1 0.186
## 7 0 1  1 0.105
## 8 1 1  1 0.047

A partir de esta conjunta podemos calcular cualquier probabilidad que nos interese. Por ejemplo, la probabilidad de que alguien compre una camisa dado que compró un pantalón es:

prob_cp <- filter(probs, p ==1 & c==1) %>% pull(prob) %>% sum
prob_cp

## [1] 0.233

prob_p <- filter(probs, p==1) %>% pull(prob) %>% sum
prob_p

## [1] 0.524

prob_cp/prob_p

## [1] 0.4446565

Como discutimos arriba, intentar estimar esta conjunta usando simples conteos de canastas no funciona, pues hay \(2^n\) posibles canastas, e incluso cuando \(n\) no es tan grande (por ejemplo 200) es un número gigantesco. Tenemos dos caminos (o una combinación de ellos): podemos hacer hipótesis acerca de esta conjunta (y checar si son apropiadas), o concentrarnos en estimar solamente algunas de sus características.

La simplificación de market basket es concentrarnos en algunas marginales que involucren a pocos artículos de esta distribución, que tienen una forma como

\[P(X_{i}=1,X_{j}=1),\] que es la probabilidad de que el conjunto \({i,j}\) aparezca en una canasta dada, o en los términos de market basket, el soporte del itemset \({i,j}\).

La búsqueda de itemsets frecuentes se traduce entonces en buscar marginales de este tipo que no involucren muchas variables y que tengan valores altos - buscamos modas en las marginales de la distribución de las canastas.

Modelo de artículos independientes

Aunque el comportamiento general de las canastas probablemente no se puede describir con modelos simples para la conjunta, puede ser útil experimentar con modelos simples para entender qué tipo de cosas podemos obervar.

En primer lugar, podemos considerar el modelo que establece que la aparición o no de cada artículo es independiente del resto:

\[P(X_1=x_1,\ldots, X_m=x_m) =\prod_m P(X_j=x_j)\]

Y adicionalmente, suponemos que la probabilidad de cada artículo es fija dada por \[P(X_j=1)=p_j\]. Entonces es fácil ver que el soporte (bajo el modelo teórico) de un conjunto de \(k\) artículos es \[P(X_{s_1}=1,X_{s_2}=1,\ldots, X_{s_k}=1)=p_{s_1}p_{s_2}\cdots p_{s_k}\]

Podemos ver qué pasa si simulamos transacciones bajo este modelo simple. Primero definimos una función para simular canastas con probabilidades dadas para los artículos

simular_transacciones <- function(nItems, nTrans, iprob){
  etiquetas <- 1:nItems
  if(!is.null(names(iprob))) {
    etiquetas <- names(iprob)
  }
  trans <- lapply(1:nTrans, function(i){
    etiquetas[which(rbinom(nItems, 1, prob = iprob) == 1)]
  })
  trans
}

Y ahora simulamos usando las proporciones que encontramos en el conjunto Groceries

probs_items <- sort(itemFrequency(Groceries))
probs_items

##                 baby food      sound storage medium 
##              0.0001016777              0.0001016777 
##     preservation products           kitchen utensil 
##              0.0002033554              0.0004067107 
##                      bags            frozen chicken 
##              0.0004067107              0.0006100661 
##            baby cosmetics            toilet cleaner 
##              0.0006100661              0.0007117438 
##            salad dressing                    whisky 
##              0.0008134215              0.0008134215 
##           make up remover                   liqueur 
##              0.0008134215              0.0009150991 
##           rubbing alcohol                hair spray 
##              0.0010167768              0.0011184545 
##             frozen fruits                     cream 
##              0.0012201322              0.0013218099 
##                     honey               decalcifier 
##              0.0015251652              0.0015251652 
##          organic products      specialty vegetables 
##              0.0016268429              0.0017285206 
##               ready soups    flower soil/fertilizer 
##              0.0018301983              0.0019318760 
##                  prosecco           organic sausage 
##              0.0020335536              0.0022369090 
##              cocoa drinks                   tidbits 
##              0.0022369090              0.0023385867 
##            pudding powder         cooking chocolate 
##              0.0023385867              0.0025419420 
##                      soap          bathroom cleaner 
##              0.0026436197              0.0027452974 
##                  cookware           potato products 
##              0.0027452974              0.0028469751 
##                      fish            snack products 
##              0.0029486528              0.0030503305 
##                 nut snack          artif. sweetener 
##              0.0031520081              0.0032536858 
##              canned fruit                     syrup 
##              0.0032536858              0.0032536858 
##               nuts/prunes          abrasive cleaner 
##              0.0033553635              0.0035587189 
##                 skin care             specialty fat 
##              0.0035587189              0.0036603965 
##                       tea                    brandy 
##              0.0038637519              0.0041687850 
##               light bulbs                   ketchup 
##              0.0041687850              0.0042704626 
##              meat spreads                       rum 
##              0.0042704626              0.0044738180 
##            male cosmetics                liver loaf 
##              0.0045754957              0.0050838841 
##               curd cheese                   cleaner 
##              0.0050838841              0.0050838841 
##                    spices                       jam 
##              0.0051855618              0.0053889171 
##                    sauces                  softener 
##              0.0054905948              0.0054905948 
##            sparkling wine                   cereals 
##              0.0055922725              0.0056939502 
##               dental care            kitchen towels 
##              0.0057956279              0.0059989832 
##  female sanitary products         finished products 
##              0.0061006609              0.0065073716 
##                   vinegar                     soups 
##              0.0065073716              0.0068124047 
##                  zwieback                   popcorn 
##              0.0069140824              0.0072191154 
##            instant coffee                      rice 
##              0.0074224708              0.0076258261 
##        liquor (appetizer)     Instant food products 
##              0.0079308592              0.0080325369 
##                    turkey    house keeping products 
##              0.0081342145              0.0083375699 
##    frozen potato products          specialty cheese 
##              0.0084392476              0.0085409253 
##                  dog food                   candles 
##              0.0085409253              0.0089476360 
##             sweet spreads     chocolate marshmallow 
##              0.0090493137              0.0090493137 
##                mayonnaise                photo/film 
##              0.0091509914              0.0092526690 
##                  pet care            condensed milk 
##              0.0094560244              0.0102694459 
##            roll products             flower (seeds) 
##              0.0102694459              0.0103711235 
##              dish cleaner            frozen dessert 
##              0.0104728012              0.0107778343 
##                      salt         canned vegetables 
##              0.0107778343              0.0107778343 
##                    liquor             spread cheese 
##              0.0110828673              0.0111845450 
##           cling film/bags               frozen fish 
##              0.0113879004              0.0116929334 
##                   mustard packaged fruit/vegetables 
##              0.0119979664              0.0130147433 
##                  cake bar         seasonal products 
##              0.0132180986              0.0142348754 
##                     pasta               canned fish 
##              0.0150482969              0.0150482969 
##                     herbs          processed cheese 
##              0.0162684291              0.0165734621 
##               soft cheese                pot plants 
##              0.0170818505              0.0172852059 
##                     flour                    dishes 
##              0.0173868836              0.0175902389 
##       semi-finished bread             baking powder 
##              0.0176919166              0.0176919166 
##        pickled vegetables                white wine 
##              0.0178952720              0.0190137265 
##            red/blush wine                 detergent 
##              0.0192170819              0.0192170819 
##               chewing gum                    grapes 
##              0.0210472801              0.0223690900 
##                  cat food             sliced cheese 
##              0.0232841891              0.0245043213 
##               hard cheese                 ice cream 
##              0.0245043213              0.0250127097 
##                      meat                       ham 
##              0.0258261312              0.0260294865 
##                 beverages             specialty bar 
##              0.0260294865              0.0273512964 
##               butter milk                       oil 
##              0.0279613625              0.0280630402 
##              frozen meals           misc. beverages 
##              0.0283680732              0.0283680732 
##                     candy       specialty chocolate 
##              0.0298932384              0.0304016268 
##                    onions          hygiene articles 
##              0.0310116929              0.0329435689 
##            hamburger meat                   berries 
##              0.0332486019              0.0332486019 
##                  UHT-milk                     sugar 
##              0.0334519573              0.0338586680 
##                   dessert  long life bakery product 
##              0.0371123538              0.0374173869 
##               salty snack                   waffles 
##              0.0378240976              0.0384341637 
##             cream cheese                white bread 
##              0.0396542959              0.0420945602 
##                   chicken         frozen vegetables 
##              0.0429079817              0.0480935435 
##                 chocolate                   napkins 
##              0.0496187087              0.0523640061 
##                      beef                      curd 
##              0.0524656838              0.0532791052 
##                    butter                      pork 
##              0.0554143366              0.0576512456 
##                    coffee                 margarine 
##              0.0580579563              0.0585663447 
##               frankfurter             domestic eggs 
##              0.0589730554              0.0634468734 
##               brown bread        whipped/sour cream 
##              0.0648703610              0.0716827656 
##     fruit/vegetable juice                 pip fruit 
##              0.0722928317              0.0756481952 
##               canned beer                newspapers 
##              0.0776817489              0.0798169802 
##              bottled beer              citrus fruit 
##              0.0805287239              0.0827656329 
##                    pastry                   sausage 
##              0.0889679715              0.0939501779 
##             shopping bags            tropical fruit 
##              0.0985256736              0.1049313676 
##           root vegetables             bottled water 
##              0.1089984748              0.1105236401 
##                    yogurt                      soda 
##              0.1395017794              0.1743772242 
##                rolls/buns          other vegetables 
##              0.1839349263              0.1934926284 
##                whole milk 
##              0.2555160142

set.seed(1299)
trans <- simular_transacciones(nItems = 169, nTrans = 10000, iprob = probs_items)
head(trans)

## [[1]]
## [1] "cocoa drinks"   "sparkling wine" "condensed milk" "chewing gum"   
## [5] "coffee"        
## 
## [[2]]
## [1] "white bread"  "citrus fruit" "sausage"      "whole milk"  
## 
## [[3]]
## [1] "packaged fruit/vegetables" "flour"                    
## [3] "white wine"                "frankfurter"              
## [5] "bottled beer"             
## 
## [[4]]
## [1] "mayonnaise" "newspapers" "soda"       "whole milk"
## 
## [[5]]
## [1] "brown bread"
## 
## [[6]]
## [1] "curd"          "bottled water" "soda"          "rolls/buns"

ap_indep <- apriori(trans, parameter = 
                      list(support=0.005, target='frequent itemsets'),
                    control = list(verbose = FALSE))
ap_indep

## set of 521 itemsets

Por ejemplo, aquí vemos algunos pares frecuentes encontrados por el algoritmo:

inspect(subset(ap_indep, support > 0.015 & size(ap_indep)==2))

##      items                                 support count
## [1]  {beef,whole milk}                     0.0151  151  
## [2]  {margarine,whole milk}                0.0160  160  
## [3]  {butter,whole milk}                   0.0159  159  
## [4]  {pork,whole milk}                     0.0172  172  
## [5]  {domestic eggs,whole milk}            0.0165  165  
## [6]  {brown bread,whole milk}              0.0186  186  
## [7]  {other vegetables,whipped/sour cream} 0.0159  159  
## [8]  {whipped/sour cream,whole milk}       0.0177  177  
## [9]  {fruit/vegetable juice,whole milk}    0.0182  182  
## [10] {pip fruit,whole milk}                0.0195  195  
## [11] {canned beer,rolls/buns}              0.0158  158  
## [12] {canned beer,other vegetables}        0.0160  160  
## [13] {canned beer,whole milk}              0.0203  203  
## [14] {newspapers,rolls/buns}               0.0150  150  
## [15] {newspapers,other vegetables}         0.0171  171  
## [16] {newspapers,whole milk}               0.0206  206  
## [17] {bottled beer,rolls/buns}             0.0151  151  
## [18] {bottled beer,other vegetables}       0.0158  158  
## [19] {bottled beer,whole milk}             0.0209  209  
## [20] {other vegetables,pastry}             0.0192  192  
## [21] {pastry,whole milk}                   0.0216  216  
## [22] {citrus fruit,soda}                   0.0151  151  
## [23] {citrus fruit,rolls/buns}             0.0167  167  
## [24] {citrus fruit,other vegetables}       0.0177  177  
## [25] {citrus fruit,whole milk}             0.0240  240  
## [26] {sausage,soda}                        0.0179  179  
## [27] {rolls/buns,sausage}                  0.0170  170  
## [28] {other vegetables,sausage}            0.0194  194  
## [29] {sausage,whole milk}                  0.0250  250  
## [30] {shopping bags,soda}                  0.0172  172  
## [31] {rolls/buns,shopping bags}            0.0150  150  
## [32] {other vegetables,shopping bags}      0.0182  182  
## [33] {shopping bags,whole milk}            0.0258  258  
## [34] {soda,tropical fruit}                 0.0180  180  
## [35] {rolls/buns,tropical fruit}           0.0193  193  
## [36] {other vegetables,tropical fruit}     0.0209  209  
## [37] {tropical fruit,whole milk}           0.0266  266  
## [38] {root vegetables,soda}                0.0184  184  
## [39] {rolls/buns,root vegetables}          0.0183  183  
## [40] {other vegetables,root vegetables}    0.0191  191  
## [41] {root vegetables,whole milk}          0.0294  294  
## [42] {bottled water,yogurt}                0.0177  177  
## [43] {bottled water,soda}                  0.0202  202  
## [44] {bottled water,rolls/buns}            0.0218  218  
## [45] {bottled water,other vegetables}      0.0227  227  
## [46] {bottled water,whole milk}            0.0287  287  
## [47] {soda,yogurt}                         0.0228  228  
## [48] {rolls/buns,yogurt}                   0.0245  245  
## [49] {other vegetables,yogurt}             0.0258  258  
## [50] {whole milk,yogurt}                   0.0334  334  
## [51] {rolls/buns,soda}                     0.0345  345  
## [52] {other vegetables,soda}               0.0344  344  
## [53] {soda,whole milk}                     0.0452  452  
## [54] {other vegetables,rolls/buns}         0.0342  342  
## [55] {rolls/buns,whole milk}               0.0462  462  
## [56] {other vegetables,whole milk}         0.0507  507

Estos pares frecuentes no se deben a asociaciones entre los artículos, sino a co-ocurrencia en las canastas. Artículos frecuentes apareceran en pares frecuentes, tríos frecuentes, etc. Comparamos por ejemplo el número de reglas encontradas para los datos reales, contra 10 simulaciones de este modelo:

sims <- lapply(1:10, function(i){
  trans <- simular_transacciones(nItems = 169, 
                                 nTrans = 10000, 
                                 iprob = probs_items)
  ap_indep <- apriori(trans, parameter = 
                        list(support=0.005, target='frequent itemsets'),
                      control = list(verbose = FALSE))
  df <- data.frame(table(size(ap_indep)))
  df$rep <- i
  df
})

sims

## [[1]]
##   Var1 Freq rep
## 1    1  116   1
## 2    2  399   1
## 3    3   13   1
## 
## [[2]]
##   Var1 Freq rep
## 1    1  118   2
## 2    2  404   2
## 3    3   15   2
## 
## [[3]]
##   Var1 Freq rep
## 1    1  122   3
## 2    2  401   3
## 3    3   13   3
## 
## [[4]]
##   Var1 Freq rep
## 1    1  119   4
## 2    2  399   4
## 3    3   15   4
## 
## [[5]]
##   Var1 Freq rep
## 1    1  120   5
## 2    2  373   5
## 3    3   19   5
## 
## [[6]]
##   Var1 Freq rep
## 1    1  116   6
## 2    2  377   6
## 3    3   17   6
## 
## [[7]]
##   Var1 Freq rep
## 1    1  117   7
## 2    2  384   7
## 3    3   14   7
## 
## [[8]]
##   Var1 Freq rep
## 1    1  118   8
## 2    2  389   8
## 3    3   19   8
## 
## [[9]]
##   Var1 Freq rep
## 1    1  118   9
## 2    2  372   9
## 3    3   17   9
## 
## [[10]]
##   Var1 Freq rep
## 1    1  117  10
## 2    2  380  10
## 3    3   15  10

Ahora vamos a comparar con el análisis de las canastas reales:

df_sims <- bind_rows(sims)
df_sims

##    Var1 Freq rep
## 1     1  116   1
## 2     2  399   1
## 3     3   13   1
## 4     1  118   2
## 5     2  404   2
## 6     3   15   2
## 7     1  122   3
## 8     2  401   3
## 9     3   13   3
## 10    1  119   4
## 11    2  399   4
## 12    3   15   4
## 13    1  120   5
## 14    2  373   5
## 15    3   19   5
## 16    1  116   6
## 17    2  377   6
## 18    3   17   6
## 19    1  117   7
## 20    2  384   7
## 21    3   14   7
## 22    1  118   8
## 23    2  389   8
## 24    3   19   8
## 25    1  118   9
## 26    2  372   9
## 27    3   17   9
## 28    1  117  10
## 29    2  380  10
## 30    3   15  10

pars <- list(supp = 0.005, target='frequent itemsets')
ap <- apriori(lista_mb, parameter = pars, control = list(verbose = FALSE))
df_obs <- data.frame(table(size(ap)))
df_obs

##   Var1 Freq
## 1    1  120
## 2    2  605
## 3    3  264
## 4    4   12

ggplot(df_sims, aes(x=Var1, y=Freq)) + geom_point() +
  geom_point(data = df_obs, colour = 'red') + labs(x = 'Tamaño')

Y vemos claramente que el modelo simple está lejos de ajustar los datos que observamos en las canastas de Groceries. Hay muchas más combinaciones frecuentes de tamaño 2 y 3 de lo que esperaríamos si los artículos se compraran independientemente, y esto indica asociaciones positivas entre artículos que nos gustaría descubrir. ¿Cómo podemos distinguir esas asociaciones?

Finalmente, comparamos los itemsets de ambos casos:

coinciden <- match(ap, ap_indep)
head(coinciden,30)
##  [1]  9 23 NA NA  1 NA  2  4 13  5 22 16 33  8  7  3 12 19 20 21 24 28 15
## [24] 11  6 40 14 30 35 27

# match produce
coinciden[500:505]

## [1] 323 324 325  NA  NA  NA

inspect(ap[500])

##     items             support    count
## [1] {pork,rolls/buns} 0.01128622 111

inspect(ap_indep[323])

##     items             support count
## [1] {pork,rolls/buns} 0.011   110

sum(!is.na(coinciden)) # contar los matches 
## [1] 491
length(ap)
## [1] 1001
length(ap_indep)
## [1] 521

Y vemos que en el análisis de datos reales estamos capturando la mayor parte de los itemsets frecuentes del modelo independiente. Estos itemsets se explican por la frecuencia simple de aparición de cada artículo.

Reglas de asociación

Aunque algunas veces lo único que nos interesa es la co-ocurrencia de artículos (por ejemplo, para entender qué artículos se podrían ver afectados por acciones en otros artículos que están en el mismo itemset frecuente), otras veces nos interesa entender qué artículos están asociados a lo largo de canastas por otros factores, como tipo de cliente, tipo de ocasión o propósito (por ejemplo, hora del día, hacer un pastel), etc.

Con este propósito podemos organizar la información de los itemsets frecuentes en términos de reglas de asociación. Un ejemplo de una regla de asociación es:

Entre las personas que compran leche y pan, un 40% compra también yogurt

Una regla de asociación es una relación de la forma \(I\to j\), donde \(I\) es un conjunto de artículos (el antecedente) y \(j\) es un artículo (el consecuente). Definimos la confianza de esta regla como

\[\hat{P}(I\to j) = \hat{P}(j|I) = \frac{n(I\cup {j})}{n(I)} = \hat{P}(I\cup j)/\hat{P}(I), \] es decir, la proporción de canastas que incluyen al itemset \(I\cup {j}\) entre las canastas que incluyen al itemset \(I\). La confianza siempre está entre 0 y 1.

Observaciones:

Por monotonicidad, si \(J\) es un conjunto de artículos más grande que \(I\) (es decir \(I\subset J\)), entonces \(n(J) \leq n(I)\): cualquier canasta que contiene a \(I\) también contiene a \(J\), y puede haber algunas canastas que contienen a \(I\) no contienen a \(J\).
Bajo el modelo de items independientes, todas las confianzas satisfacen \(P(I\to j)=P(j)\) (la confianza simplemente es la probabilidad de observar el artículo \(j\), independientemente del antecedente).
Confianza alta no necesariamente significa una asociación de los items: si el consecuente \(j\) tiene soporte alto, entonces podemos obtener confianza alta aunque no haya asociación.

Ejemplo

En nuestro ejemplo anterior, el soporte de {whole milk,yogurt} es de 0.0560, el soporte de {whole milk} es 0.2555, así que la confianza de la regla \(whole milk \to yogurt\) es \(\frac{0.0560}{0.2555}=\) 0.22

Podemos usar la confianza para filtrar reglas que tienen alta probabilidad de cumplirse:

Ejemplo

pars <- list(supp = 0.01, confidence = 0.20, target='rules', 
             ext = TRUE, minlen = 2)
reglas <- apriori(lista_mb, parameter = pars)

Podemos examinar algunas de las reglas:

inspect(reglas[1:10,])

##      lhs                rhs                support    confidence
## [1]  {hard cheese}   => {whole milk}       0.01006609 0.4107884 
## [2]  {butter milk}   => {other vegetables} 0.01037112 0.3709091 
## [3]  {butter milk}   => {whole milk}       0.01159126 0.4145455 
## [4]  {ham}           => {whole milk}       0.01148958 0.4414062 
## [5]  {sliced cheese} => {whole milk}       0.01077783 0.4398340 
## [6]  {oil}           => {whole milk}       0.01128622 0.4021739 
## [7]  {onions}        => {other vegetables} 0.01423488 0.4590164 
## [8]  {onions}        => {whole milk}       0.01209964 0.3901639 
## [9]  {berries}       => {yogurt}           0.01057448 0.3180428 
## [10] {berries}       => {other vegetables} 0.01026945 0.3088685 
##      lhs.support lift     count
## [1]  0.02450432  1.607682  99  
## [2]  0.02796136  1.916916 102  
## [3]  0.02796136  1.622385 114  
## [4]  0.02602949  1.727509 113  
## [5]  0.02450432  1.721356 106  
## [6]  0.02806304  1.573968 111  
## [7]  0.03101169  2.372268 140  
## [8]  0.03101169  1.526965 119  
## [9]  0.03324860  2.279848 104  
## [10] 0.03324860  1.596280 101

En la siguiente tabla, lhs.support es el soporte del antecedente (lhs = left hand side). Agregamos también el error estándar de la estimación de confidence (que es una proporción basada en el número de veces que se observa el antecedente):

df_1 <- sort(reglas, by = 'confidence') %>%
  DATAFRAME 
df_2 <- df_1 %>% select(LHS, RHS, lhs.support, confidence, support) %>%
  head(100) %>%
  mutate(lhs.base = num_canastas*lhs.support) %>%
  mutate(conf.ee = sqrt(confidence*(1-confidence)/lhs.base)) %>%
  mutate_if(is.numeric, funs(round(., 2))) 
DT::datatable(df_2)

Observaciones:

Nota que estas tres cantidades están ligadas en cada canasta por \(lhs.support\times confidence = support\). Usa un argumento de probabilidad condicional para mostrarlo.
Muchas de las reglas con confianza alta tienen como consecuente un artículo de soporte alto (por ejemplo, whole milk), como explicamos arriba. Nótese también que las reglas con confianza más alta tienden a tener soporte bajo. Esto lo discutiremos más adelante.

Ejercicio

Para un mismo consecuente (por ejemplo whole milk), examina cómo varían los valores de confidence. ¿A qué crees que se deba esto?

Es natural que artículos frecuentes ocurran en muchas canastas juntas, es decir, que reglas formadas con ellas tengan confianza relativamente alta. Por ejemplo, la regla pan -> verduras podría tener confianza y soporte alto, pero esto no indica ninguna asociación especial entre estos artículos. La razón puede ser que verduras es un artículo muy común

Podemos refinar las reglas de asociación considerando qué tan diferente es \(P(j|I)\) de \(P(j)\). La primera cantidad es la probabilidad de observar el item \(j\) bajo la información de que la canasta contiene a \(I\). Si esta cantidad no es muy diferente a \(P(j)\), entonces consideramos que esta regla no tiene mucho interés.

El lift o intéres de una regla \(I\to j\) se define como \[L(I\to j) = \frac{\hat{P}({j}|I)}{\hat{P}({j})},\] es decir, la confianza de la regla \(I\to j\) dividida entre la proporción de canastas que contienen \(j\).

En nuestro ejemplo, veamos dos reglas con interés muy distinto:

df_1 <- arrange(df_1, desc(lift))
df_1[c(4, nrow(df_1)),] %>% select(LHS, RHS, lhs.support, confidence, lift)

##        LHS               RHS lhs.support confidence      lift
## 4   {beef} {root vegetables}  0.05246568  0.3313953 3.0403668
## 231 {soda}      {whole milk}  0.17437722  0.2297376 0.8991124

La primera regla tiene un interés mucho más alto que la segunda, lo que indica una asociación más importante entre los dos artículos.

Observaciones

Cuando decimos que un grupo de artículos están asociados, generalmente estamos indicando que forma alguna regla de asociación con lift alto.
En principio también podría haber reglas con lift muy por debajo de uno, y eso también indica una asociación (por ejemplo coca y pepsi). Pero el método de itemsets frecuentes no es muy apropiado para buscar estas reglas, pues precisamente esas reglas tienden a tener soporte y confianza bajas.
El valor del lift también puede escribirse (deméstralo) como \[ \frac{\hat{P}(I\cup\{j\})}{\hat{P}(I)\hat{P}({j})},\] Cuando los artículos son independientes, esta cantidad está cercana a 1. Es una medida de qué tan lejos de independencia están la ocurrencia de los itemsets \(I\) y \(j\).

Ejemplo

df_1 <- sort(reglas, by = 'lift') %>%
  DATAFRAME
df_1

En esta tabla, lhs.support es el soporte del antecedente (lhs = left hand side):

df_2 <- df_1 %>% select(LHS, RHS, lhs.support, lift, confidence, support) %>%
  head(100) %>%
  mutate_if(is.numeric, funs(round(., 2))) 
DT::datatable(df_2)

Las reglas de asociación se calculan comenzando por calcular los itemsets frecuentes según el algoritmo a priori explicado arriba. Para encontrar las reglas de asociación hacemos:

Para cada itemset frecuente \(f\), tomamos como candidatos a consecuentes los artículos \(i\) de \(f\)
Si la confianza \(\frac{\hat{P}(I)}{\hat{P}(I-\{j\})}\) es mayor que la confianza mínima, agregamos la regla de asociación \(I\to j\).

Con este proceso encontramos todas las reglas \(I\to j\) tales que \(I\cup\{j\}\) es un itemset frecuente.

Dificultades en el análisis de canastas

El análisis de canastas es un método rápido y simple que nos da varias maneras de explorar las relaciones entre los artículos. Sin embargo, hay varias dificultades en su aplicación.

Número de reglas y itemsets

Muchas veces encontramos un número muy grande de itemsets o reglas. Hay varias maneras de filtrar estas reglas según el propósito. Si filtramos mucho, perdemos reglas que pueden ser interesantes. Si filtramos poco, es difícil entender globalmente los resultados del análisis.

Un punto de vista es producir una cantidad de reglas para procesar posteriormente con queries: por ejemplo, si nos interesa entender las relaciones de berries con otros artículos, podemos filtrar las reglas encontradas y examinarlas más fácilmente.

Cortes estrictos en el filtrado

Cuando seleccionamos valores mínimos de soporte, confianza y/o lift, estas decisiones son más o menos arbitrarias. Distintos analistas pueden llegar a resultados distintos, incluso cuando el propósito del análisis sea similar, y en ocasiones hay que iterar el análisis para encontrar valores adecuados que den conjuntos razonables con resultados interesantes. Este último concepto es subjetivo.

Redundancia de reglas

Existe alguna redundancia en las reglas que encontramos. Por ejemplo, podríamos tener {yogurt, berries} -> {whipped cream}, pero también {yogurt} -> {whipped cream}. Este traslape de reglas hace también difícil entender conjuntos grandes de reglas.

Variabilidad de medidas de calidad

Un problema del análisis clásico de soporte-confianza-lift es la variabilidad de las estimaciones de confianza y lift.

Cuando comenzamos poniendo valores de soporte y confianza relativamente bajos, encontramos muchas reglas
Intentamos muchas veces filtrar u ordenar por lift, para considerar las reglas más interesantes
Sin embargo, encontramos entonces que muchas reglas de lift y/o confianza altas son aquellas que tienen soporte bajo y consecuentes poco frecuentes. Como veremos más adelante, esto se debe muchas veces a error de estimación. Los valores más grandes de lift generalmente son sobreestimaciones, por la naturaleza del análisis basado en cortes.
Si regresamos a incrementar soporte y confianza, potencialmente perdemos reglas interesantes.

Podemos usar el modelo de independencia Veamos cómo se comportan confianza y lift para el modelo donde no hay asociaciones. Obsérvese que en este modelo todas las confianzas son iguales a la frecuencia del consecuente, y todos los valores de lift son 1:

sims_reglas <- lapply(1:10, function(i){
  trans <- simular_transacciones(nItems = 169, 
                                 nTrans = 10000, 
                                 iprob = probs_items)
  ap_random <- apriori(trans, parameter = 
                         list(support=0.002, confidence = 0.0, 
                              target='rules', ext = TRUE, minlen = 2),
                       control = list(verbose = FALSE))
  ap_random
})

Y notamos que conforme el soporte de la regla es más bajo, hay más variabilidad en las estimaciones del confianza y lift. En este caso utilizamos

plotly_arules(subset(sims_reglas[[4]], rhs %pin% 'whole milk'), 
              measure=c('support','confidence'), shading = 'lift')

plotly_arules(subset(sims_reglas[[4]], rhs %pin% 'whole milk'), 
              measure=c('support','lift'), shading = 'confidence')

El valor de confianza y de lift puede ser altamente variable para reglas con soporte bajo. Podemos tomar dos caminos:

Cuando hagamos el soporte más bajo, incrementamos el valor de lift mínimo. Esto evita que obtengamos demasiadas reglas que no representan interacciones reales entre los artículos.
Podemos usar otras medidas que tomen en cuenta la variabilidad de las estimaciones. Por ejemplo, hyper-lift y hyper-confidence están basados en modelos simples (como el que vimos arriba), que filtran aquellos valores de calidad que están en las colas de las distribuciones de los modelos simples.

Otras medidas de calidad de reglas

Hay una gran cantidad de medidas de interés de reglas que se han propuesto desde que se usa el análisis de canasta. Aquí discutimos hyper-lift y hyper-confidence, que toman en cuenta el soporte de las reglas para proponer puntos de corte [@hyper].

Explicamos aquí el hyper-lift para una regla \(i\to j\).

Consideramos el modelo de independencia (lo pensamos como el modelo nulo), fijando las probabilidades de ocurrencia de los artículos según los datos (como hicimos en los ejemplos de arriba) y el número de transacciones. Bajo este modelo, el número de ocurrencias \(X_{\{i,j\}}\) de el itemset \(\{i,j\}\) es una variable aleatoria con distribución conocida (binomial). Esta distribución representa la variación que podemos observar en los conteos de \(\{i,j\}\) bajo distintas muestras de transacciones del mismo tamaño.

La idea básica del hyperlift es comparar el conteo \(n(\{i,j\})\) con la cola superior de la distribución de \(X_{i,j}\) bajo el supuesto de independencia, poniendo

\[HL(I\to j) = \frac{n(\{i,j\})}{Q_\delta (X_{i,j})},\]

donde \(Q_\delta\) es tal que \(P(X_{i,j} > Q_\delta (X_{i,j}))\approx \delta\). Tomamos por ejemplo \(\delta=0.99\). De esta forma, \(HL>1\) sólo cuando el conteo observado \(n(\{i,j\})\) está en la cola superior del conteo bajo la hipótesis nula. Esto toma en cuenta la variabilidad de los conteos (que es grande en términos relativos cuando el soporte es bajo).

Observaciones:

El modelo de independencia que se usa en el paquete arules es una variación del que vimos aquí, ver los detalles en [@hyper].
Los valores de hyper-lift no son realmente comparables a los de lift, son dos medidas de calidad de asociación diferentes, pero similares en cuanto a lo que quieren capturar.

trans <- simular_transacciones(nItems = 169, 
                               nTrans = 10000, 
                               iprob = probs_items)
ap_random <- apriori(trans, parameter = 
                       list(support=0.001, confidence = 0.10, 
                            target='rules', ext = TRUE, minlen = 2),
                     control = list(verbose = FALSE))
agregar_hyperlift <- function(reglas, trans){
  quality(reglas) <- cbind(quality(reglas), 
                           hyper_lift = interestMeasure(reglas, measure = "hyperLift", 
                                                        transactions = trans))
  reglas
}
ap_random <- agregar_hyperlift(ap_random, trans)

Podemos ver cómo se relacionan lift y hyper_lift. Vemos claramente que la gran mayoría de reglas obtenidas ahora tienen hyper-lift menor que uno

plot(ap_random, measure=c('lift','hyper_lift'), shading = 'support')

Cortando en un valor relativamente bajo, vemos que nos deshacemos correctamente de casi todas las reglas:

length(ap_random)
## [1] 3625
length(subset(ap_random, lift > 1))
## [1] 2742
length(subset(ap_random, hyper_lift > 1))
## [1] 127

Ahora aplicamos a los datos reales

pars <- list(supp = 0.002, confidence = 0.10, target='rules', 
             ext = TRUE, minlen = 2)
reglas <- apriori(lista_mb, parameter = pars,
                  control = list(verbose=FALSE))
reglas <- agregar_hyperlift(reglas, Groceries)
length(reglas)
## [1] 8332

Vemos que podemos cortar niveles de hyper-lift donde obtenemos reglas de soporte relativamente alto.

plot(reglas, measure=c('lift','hyper_lift'), shading = 'support')

Si cortamos en valores que dan un número similar de reglas, por ejemplo:

filtradas_hl <- subset(reglas, hyper_lift > 2)
filtradas_lift <- subset(reglas, lift > 3.7)
length(filtradas_hl)

## [1] 445

length(filtradas_lift)

## [1] 480

Vemos que las reglas cortadas con hyperlift tienen mejores valores de soporte

qqplot(quality(filtradas_hl)$support, 
       quality(filtradas_lift)$support, xlab = 'Hyper-lift',
       ylab = 'Lift', main = 'Soporte')

Sin embargo, la distribución de valores de lift no es tan diferente:

qqplot(quality(filtradas_hl)$lift, 
       quality(filtradas_lift)$lift,
       xlab = 'Hyper-lift',
       ylab = 'Lift', main = 'Lift')

En resumen, al descartar aquellas reglas con soporte bajo que tienen lift alto por azar, encontramos reglas de mejor calidad en general.

Selección de reglas

Ahora discutiremos cómo seleccionar itemsets frecuentes y reglas.

Filtrar con todos estos criterios (soporte, confianza, soporte del antecedente, lift) no es simple, y depende de los objetivos del análisis. Recordemos también que estos análisis están basados justamente en cortes “duros” de los datos, más o menos arbitrarios, y por lo tanto pueden los resultados son variables.

Hay varias maneras de conducir el análisis. Dos tipos útiles son:

Itemsets de alta frecuencia: en este enfoque buscamos reglas con soporte y confianza relativamente altos. Generalmente están asociados a productos muy frecuentes, y nos indica potencial de interacción entre los artículos. Este análisis es más una reexpresión de la información contenida en los itemsets frecuentes. En este caso, podemos filtrar con soporte alto, para evitar estimaciones ruidosas (por ejemplo, soporte mínimo de 300 canastas).
Interacciones altas: en este enfoque donde buscamos entender nichos. Consideramos valores de soporte y confianza más bajos, pero con valores de lift alto. Este análisis es más útil para entender, por ejemplo, propósitos de compras, convivencia de artículos, tipos de comprador, etc.
Colección de reglas para hacer querys: la colección de reglas puede ser más grande, e incluir por ejemplo resultados de distintas corridas de market basket (incluyendo los dos enfoques de arriba). Las reglas se examinan seleccionando antecedentes o consecuentes, valores altos de soporte, etc, según la pregunta particular que se quiera explorar.

Ejemplo: canastas grandes

Para entender las canastas grandes, podemos variar valores de soporte y confianza para encontrar un número manejable de reglas.

pars <- list(support = 0.02,
             confidence = 0.20,
             minlen = 2,
             target='rules', ext = TRUE)
reglas_1 <- apriori(lista_mb, parameter = pars)

Esta elección de parámetros resulta en 72. Podemos ordenar por hyper-lift:

plotly_arules(reglas_1, colors=c('red','gray'))

reglas_1 <- agregar_hyperlift(reglas_1, lista_mb)
DT::datatable(DATAFRAME(sort(reglas_1, by='hyper_lift')) %>%
                mutate_if(is.numeric, funs(round(., 3))))

Observaciones: conforme bajamos en esta tabla (ordenada por soporte), las estimaciones de confianza y lift son menos precisas.

Búsqueda de reglas especializadas

Otra manera de usar este análisis es intenando buscar asociaciones más fuertes (lift o hyper-lift más alto), aún cuando sacrificamos soporte. Por su naturaleza, este tipo de análisis puede resultar en reglas más ruidosas (malas estimaciones de confianza y lift), pero es posible filtrar valores más altos de estas cantidades para encontrar reglas útiles.

Comenzamos con un soporte y confianza más bajas

pars <- list(support = 0.001,
             confidence = 0.1,
             minlen = 2,
             target='rules', ext = TRUE)
b_reglas <- apriori(lista_mb, parameter = pars)
b_reglas <- agregar_hyperlift(b_reglas, lista_mb)

Y ahora filtramos con valores más grandes de hyper-lift. Podemos filtrar adicionalmente con lhs.support para obtener reglas que aplican con más frecuencia:

b_reglas

## set of 32783 rules

b_reglas_lift <- subset(b_reglas, 
                        hyper_lift > 2.5 & size(b_reglas) < 4 &
                          lhs.support > 0.01)
b_reglas_lift <- sort(b_reglas_lift, by = 'hyper_lift')
DT::datatable(DATAFRAME(b_reglas_lift)  %>%
                mutate_if(is.numeric, funs(round(., 3))))

Visualización de asociaciones

Tener una visión amplia del market basket analysis es difícil (típicamente, funciona mejor como un resultado al que se le hacen querys, o uno donde filtramos cuidadosamente algunas reglas que puedan ser útiles). Así que muchas veces ayuda visualizar los pares con asociación alta:

Construimos todas las reglas con un antecedente y un consecuente.
Filtramos las reglas con hyper-lift relativamente alto (por ejemplo > 1.5, pero hay que experimentar).
Representamos como una gráfica donde los nodos son artículos, y las aristas son relaciones de lift alto.
Usamos algún algoritmo para representar gráficas basado en fuerza, usando como peso el lift.

Ejemplo

En nuestro caso, podríamos tomar (ajustando parámetros para no obtener demasiadas reglas o demasiado pocas)

b_reglas_lift <- subset(b_reglas, 
                        hyper_lift > 1.75 & confidence > 0.1)
reglas_f <- subset(b_reglas_lift, size(b_reglas_lift)==2)

library(tidygraph)

## 
## Attaching package: 'tidygraph'

## The following object is masked from 'package:stats':
## 
##     filter

library(ggraph)
df_reglas <- reglas_f %>% DATAFRAME %>% rename(from=LHS, to=RHS) %>% as_data_frame
df_reglas$weight <- log(df_reglas$hyper_lift)
graph_1 <- as_tbl_graph(df_reglas) %>%
  mutate(centrality = centrality_degree(mode = "all")) 

ggraph(graph_1, layout = 'fr', start.temp=100) +
  geom_edge_link(aes(alpha=lift), 
                 colour = 'red',
                 arrow = arrow(length = unit(4, 'mm'))) + 
  geom_node_point(aes(size = centrality, colour = centrality)) + 
  geom_node_text(aes(label = name), size=4,
                 colour = 'gray20', repel=TRUE) +
  theme_graph()

Para gráficas más grandes, es mejor usar software especializado para investigar las redes que obtenemos (como Gephi):

library(readr)
b_reglas_lift <- subset(b_reglas, 
                        hyper_lift > 1.75 & confidence > 0.05)
reglas_f <- subset(b_reglas_lift, size(b_reglas_lift)==2)
# write_csv(reglas_f %>% DATAFRAME %>% rename(source=LHS, target=RHS) %>%
#             select(-count),
#           path='./salidas/reglas.csv')

Para el análisis de canastas grandes:

reglas_f2 <- subset(reglas_1, hyper_lift > 1.3, confidence > 0.4)
df_reglas <- reglas_f2 %>% DATAFRAME %>% rename(from=LHS, to=RHS) %>% as_data_frame
df_reglas$weight <- log(df_reglas$hyper_lift)
graph_1 <- as_tbl_graph(df_reglas) %>%
  mutate(centrality = centrality_degree(mode = 'all'))

ggraph(graph_1, layout = 'fr', start.temp=100) +
  geom_edge_link(aes(alpha=hyper_lift), colour = 'red',arrow = arrow(length = unit(4, 'mm'))) + 
  geom_node_point(aes(size = centrality, colour = centrality)) + 
  geom_node_text(aes(label = name), size=4,
                 colour = 'gray20', repel=TRUE) +
  theme_graph()

Otras aplicaciones

Análisis de tablas de variables categóricas: podemos considerar una tabla con varias variables categóricas. Una canasta son los valores que toman las variables. Por ejemplo, podríamos encontrar reglas como {hogar = propio, ocupación=profesional} -> ingreso = alto. Puedes ver más de este análisis en [@esl], por ejemplo, sección 14.2.
Conceptos relacionados: si los artículos son palabras y las canastas documentos (tweets, por ejemplo), este tipo de análisis (una vez que eliminamos las palabras más frecuentes, que no tienen significado como artículos, preposiciones, etc.), puede mostrar palabras que co-ocurren para formar conceptos relacionados (por ejemplo, )
Plagiarismo: si los artículos son documentos y los canastas oraciones, el análisis de canastas puede encontrar documentos que contienen las mismas oraciones. Si varias canastas (oraciones) “contienen” los mismos artículos (documentos), entonces esas oraciones son indicadores de plagio

Tarea

Considera los datos de datos/recetas. Lee los datos, asegúrate que puedes filtrar por tipo de cocina, y que puedes aplicarles la función apriori de arules (o cualquier otra herramienta que estés utilizando). Calcula la frecuencia de todos los artículos (ingredientes). El resto de este ejercicio lo haremos a principio de la siguiente clase. Acerca de los datos: Cada receta es una canasta, y los artículos son los ingredientes que usan. Puedes consultar el artículo original aquí.
Haz algunos experimentos el ejemplo @ref(ejemplo-canastas) que vimos en clase: incrementa/decrementa hyperlift, incrementa/decrementa soporte. ¿Qué pasa con las gráficas resultantes y el número de reglas?
(Opcional) Muchas veces el análisis de canastas puede hacerse con una muestra de transacciones. Leer secciones 6.4.1 a 6.4.4 de [@mmd].